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HOW  ACCURATE  IS  MATHEMATICAL  FIRE  MODELING? 


Henri  E.  Mitler  and  John  A.  Rockett 
Abstract 

It  is  important  to  be  able  to  predict  the  development  of  a fire  in  an 
enclosure  of  arbitrary  complexity.  A mathematical  model  valid  for  a single 
room,  with  multiple  vents  and  objects  in  it  has  been  developed.  The  fifth 
version  of  the  model  has  just  been  completed  ; it  is  the  Harvard  Computer  Fire 
Code  V,  or  Mark  5 for  short.  In  1977,  Factory  Mutual  Research  Corporation 
carried  out  a series  of  eight  well-instrumented  full  scale  room  fires,  against 
which  the  single  room  model  can  be  tested.  The  test  fire  room  was  2.4  m x 
9.6  m x 2.4  m high,  with  an  open  doorway;  a slab  of  polyurethane  foam  in  one 
corner,  and  a polyurethane  foam  target  in  a facing  corner.  The  primary  slab 
was  ignited  at  its  center,  and  the  fire  followed.  The  other  tests  were 
variants  of  this  one.  We  compare  the  results  of  the  calculations  with  the 
results  of  two  of  the  experiments:  the  standard  configuration  and  the  case 
with  a window  replacing  the  doorway.  The  model  "predictions"  are  in  good  to 
excellent  agreement  for  most  of  the  variables.  The  disagreements  are 
discussed,  and  it  is  found  that  the  most  probable  causes  are:  failure  to  take 
heating  of  the  bottom  gas  layer  into  account,  inadequacy  of  the  burnout 
algorithm,  and  the  lack  of  understanding  of  the  CO-production  mechanism. 


ifc 

This  article  is  based  on  a paper  presented  at  the  US-USSR  Scientific  Exchange 
Meeting,  Tbilisi,  USSR,  July  1981,  which  has  appeared  in  Russian.  The  long 
delay  in  publication  is  partly  due  to  an  agreement  with  the  USSR  to  wait  for 
preliminary  Russian  publication.  There  has  since  then  been  substantial 
development  of  Mark  5 (see  ref.  [21]). 
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room  fires. 


1 . INTRODUCTION 

The  growth  of  a fire  in  a building  is  an  extremely  complex  phenomenon. 
Nevertheless,  with  modern  scientific  tools,  it  should  be  amenable  to  analytic 
prediction.  Perhaps  the  first  research  in  predicting  fire  behavior  was  that 
of  Ingberg  [1]  at  the  National  Bureau  of  Standards.  He  considered  the  effect 
of  amount  of  fuel  (fuel  loading)  on  fully-involved  room  fires;  Kawagoe  [2] 
improved  the  analysis  by  taking  ventilation  into  account.  The  last  decade  has 
seen  a rapid  development  in  the  understanding  of  fire  physics  and  in  synthe- 
sizing and  quantifying  this  understanding  in  analytic  models.  These  advances 
will  not  be  reviewed  here.  Good  bibliographies  can  be  found  in  Mitler  [3]  and 
Babrauskas  and  Williamson  [4],  Much  has  also  been  learned  about  physically 
modeling  fire  (Heskestad  [5]).  However,  this  discussion  will  be  confined  to 
mathematical  modeling. 

A significant  feature  of  recent  fire  models  is  their  attention  to  the 
initial  growth  of  fire  from  a small  ignition  source  to  full  room  involvement 
(flashover),  and  beyond.  This  is  in  response  to  the  belief  of  many  fire 
protection  engineers  that  fire  can  be  controlled  before  flashover.  Emphasis 
on  the  early  fire  has  been  possible  because  the  physics  of  fire  is  now  much 
better  understood  than  it  was  only  ten  or  fifteen  years  ago. 

The  pre-f  lashover  fire  is  inherently  more  complicated  than  the  post- 
flashover  fire,  but  the  development  of  modern  computers  permits  the  solution 
of  the  equations,  at  least  in  some  restricted  form. 
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There  are  basically  two  kinds  of  mathematical  models.  First,  "field"  (or 
continuum)  models,  where  the  Navier-Stokes  equations  are  used  for  the  fluid 
flows,  together  with  equations  describing  the  combustion  and  radiation. 
Second,  "zone"  models,  which  employ  much-simplified  equations  to  describe  each 
of  the  several  phenomena  which  can  be  conceptually  isolated  in  a zone.  The 
latter  seem  (for  the  moment)  more  practical,  since  the  solutions  of  their 
equations,  even  in  the  time-dependent  case,  are  literally  orders  of  magnitude 
faster  than  the  solution  of  the  more  nearly  exact  field  equations.  Zone 
models  are  able  to  represent  a much  more  varied  physical  content  than  field 
models  but  at  the  cost  of  much  less  spatial  resolution.  There  is  a degree  of 
artistry  (or  arbitrariness)  involved  in  zone  "modeling",  since  the  models 
representing  each  phenomenon  must  be  simultaneously  reasonably  accurate  yet 
simple  enough  to  be  described  by  equations  which  are  readily  solved. 

Numerous  studies  have  been  carried  out  with  zone  models;  of  particular 
note  is  the  work  of  Quintiere  [6],  MacArthur  [7],  Pape  [8],  and  Emmons, 
Mitler,  and  Trefethen  [9].  A survey  of  recent  developments  in  modeling  is 

j||f 

given  by  Friedman  [10]. 

Ultimately,  it  should  be  possible  to  predict  the  development  of  a fire  in 
an  enclosure  of  arbitrary  complexity.  What  has  been  achieved  so  far  is  a fair 
simulation  of  what  happens  in  a single  room  with  maltiple  vents  and  one  or 
more  objects.  Several  multi-room  models  are  also  available  of  which  the  most 
elaborate  is  probably  that  of  Tanaka  [11].  His  model,  however,  while  allowing 
for  multiple  rooms  at  multiple  levels,  assumes  a fire  with  a prescribed  energy 
release  rate  as  a function  of  time.  There  are  strong  similarities  among  the 

Also  see  refs.  [22]  and  [23] . 
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single  room  zone  models  currently  in  use  although  they  differ  in  detail.  In 
this  paper  attention  is  focused  on  the  simulation  developed  at  Harvard  over 
the  past  seven  years.  The  latest  major  revision  of  this  computer  program  is 
called  the  Harvard  Fire  Code  V,  or  Mark  5 for  short.  It  simulates  fire  in  a 
single  room  vented  directly  to  ambient  air.  It  allows  for  up  to  five  vents 
and  five  combustible  objects.  Each  object  can  heat  up  and  possibly  ignite. 
The  simulation  allows  for  one  of  three  possible  initial  fires:  a burner 

fire — i.e.,  one  of  constant  radius  and  heat  release  rate,  a pool  fire — i.e., 
one  of  constant  radius,  and  a growing  fire — one  ignited  over  a small  area, 
which  then  spreads. 

This  paper  will  describe  the  principal  features  of  the  Harvard  Mark  5 
simulation  and  then  compare  its  predictions  for  two  fires  with  the 
experimentally  observed  behavior. 

2.  SCOPE  OF  THE  MODEL 


We  will  describe  some  of  the  pheonomena  which  take  place  when  a fire  is 
ignited  in  a room,  below.  This  has  been  done  in  more  detail  elsewhere  (see, 
for  example,  ref.  [3]).  Most  of  the  important  physical  effects  involved  in  a 
room  fire  have  been  explicitly  taken  into  account  in  Mark  5,  to  one  or  another 
degree  of  approximation  [3,9-13] . The  phenomena  not  included  are  mainly 
associated  with  post-f lashover  behavior  or  events  immediately  preceding 
flashover.  Although  classed  as  a zone  model,  the  zones  of  Mark  5 are  not 
necessarily  clearly  delineated  subdivisions  of  the  room  volume . They  are: 
(1)  the  objects,  (2)  the  fire  plumes  over  any  objects  which  are  burning, 
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(3)  the  layer  of  heated  gas  trapped  under  the  ceiling  (assumed  to  have  a 
sharply  delineated  lower  surface  although,  in  practice,  this  may  not  be  the 
case),  (4)  the  ceiling  and  that  portion  of  the  walls  in  contact  with  hot  gas 
layer,  (5)  the  floor  and  lower  walls,  (6)  the  vents.  Note  that  the  vents  are 
a phenomenological  zone  but  physically  a boundary  between  the  inside  and 
outside  of  the  room.  Some  of  the  phenomena  associated  with  the  zones  are: 

(a)  entrainment  of  air  into  the  hot,  buoyant  plume  produced  by  the  fire, 

(b)  the  production  of  soot  and  gaseous  combustion  products  and  depletion  of 
oxygen  in  the  plumes,  (c)  the  (transient)  formation  of  the  hot  layer  under  the 
ceiling,  (d)  radiative  energy  exchange  between  the  fire  plume,  hot  layer, 
upper  and  lower  walls,  ceiling  and  floor  and  with  the  space  outside  the  vents, 

(e)  the  heating  of  the  ceiling  and  upper  walls  and  the  objects  in  the  room, 

(f)  ignition  of  objects  not  initially  burning,  either  through  auto-ignition  or 

by  flame  contact,  (g)  mass  flow  through  the  vents.  Some  situations  not 

currently  included  are:  (1)  vents  located  in  the  ceiling,  (2)  external  winds, 

(3)  mixing  at  the  vents  between  the  hot  exiting  fire  gases  and  entering 

ambient  air,  together  with  effects  associated  with  this  mixing  on  the  lower, 
cool  gas  layer,  (4)  the  effect  of  pre-ignition  pyrolysis,  burning  in  vitiated 
air  and  upper  layer  flashover,  (3)  the  effect  of  sprinklers  or  other  extin- 
guishers. Although  CO  production  is  calculated,  the  present  method  for  doing 
this  reflects  our  very  imperfect  understanding  of  the  production  mechanism  and 
is  based  on  very  limited  data  correlations.  Finally,  the  effects  of  initial 
stable  gas  stratification,  gas  compression  and  finite  transit  times  have  been 
assumed  to  be  of  minor  importance.  This  last  could  become  a problem  were  Mark 
5 used  for  estimating  the  performance  of  fire  detectors  in  high  ceilinged 
rooms . 
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In  the  following  paragraphs,  a brief  outline  of  the  equations  defining 
the  model  will  be  given.  This  has  already  been  done  in  some  detail  by  Emmons 
[12],  and  in  more  detail  in  reference  [9] . The  model  has,  of  course,  evolved 
since  these  reports  were  written,  and  the  more  comprehensive  set  of  equations 
is  found  in  references  [3]  and  [13] . 

Figure  1 shows  the  situation  to  be  simulated  a few  mintues  into  the  fire, 
and  some  of  the  room  dimensions  (h)  and  selected  variables,  such  as  the  mass 
flows  (m)  and  radiative  fluxes  (<j>). 

3.  MODELS  OF  THE  PHYSICAL  PROCESSES 


3.1  The  Burning  Object 

As  stated  above,  Mark  5 allows  for  three  types  of  fires:  a gas  burner,  a 

pool  fire  or  a growing  fire.  The  growing  fire  algorithm  is  the  most  general, 
the  other  two  are  specializations  of  this.  The  growing  fire  model  is  for  a 
centrally  ignited  horizontal  slab  of  material.  It  was  experimentally  observed 
that,  for  foamed  polyurethane,  the  fire  radius  increases  exponentially  with 
time  [14].  We  assume  that  the  spread  rate  depends  on  the  net  impinging  flux 
in  such  a way  that  the  experimentally-observed  exponential  growth  rate  for  a 
small  fire  is  reproduced.  Reradiation  and  convective  heating  of  the  burning 
surface  are  taken  into  account  in  calculating  the  net  impinging  heat  flux. 
This  leads  to  the  spread-rate  expression  for  the  fire  radius, 

r r2 

R = AC  (1  + T + o ) , R<  °*95  R ( 1 ) 

^ 2 3 ' - m 
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where 


C = <}>/aT^ 


(2) 


and  A is  a spread-rate  parameter  specific  to  the  material.  <j>  is  the  net  heat 
flux  impinging  on  the  pyrolyzing  surface,  and  is  the  flame  temperature.  It 
will  be  seen  that  the  right  hand  side  of  equation  (1)  represents  the  first 
three  terms  in  the  expansion  of  log(l-C).  This  was  done  in  order  to  avoid  the 
singularity  at  C = 1. 

When  the  radius  R approaches  the  input  maximum  radius,  Rm,  the  spread 
rate  evidently  must  go  down,  and  so  we  take 

R = (R  - Rl/10,  R > 0.95  R (3) 

m 1 m 

R is  integrated  to  find  the  radius,  R(t). 

The  pyrolysis  rate  is  then  given  by 

raf  = ~ tiR2(J./H  (A) 

f Y v 

where  Hy  is  the  heat  of  pyrolysis  (or  vaporization).  The  negative  sign  is 
included  because  m^(t)  is  the  mass  of  the  fuel. 

An  important  point  is  that  <j>  includes  all  the  radiation  incident  on  the 
surface,  including  that  from  the  hot  layer,  walls,  and  its  own  flame  and  other 
fires.  This  feedback  produces  an  acceleration  of  burning  in  the  room, 

relative  to  burning  in  the  open  where  only  radiation  from  the  objects  own 
flame  would  be  present. 
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Combustion  in  the  diffusion  flames  found  in  fires  is  typically  quite 
inefficient.  Only  a fraction  of  the  pyrolyzed  fuel  in  fact  burns,  the  rest 
going  off  as  soot,  CO,  and  unburned  hydrocarbons.  Furthermore,  when  the  hot 
upper  gas  layer  is  so  low  that  it  keeps  the  fire  plume  from  entraining  much 
air,  the  fire  becomes  starved  for  oxygen,  and  the  burning  rate  may  go  lower 
still.  Thus  the  burning  rate  is 

= min  (-  x mf  , rWy)  (5) 

where  m^  is  the  rate  of  entrainment  of  air  into  the  plume,  y is  the  effective 
air/fuel  mass  ratio  in  burning,  and  x is  the  fraction  of  pyrolyzed  fuel  which 
actually  burns  (completely).  This  is  found  from  experiment  to  be  about 
14.45.  The  chemical  energy  released  then  in  either  case  is 

q = mbHc>  (6) 

where  H = the  heat  of  combustion 
c 

The  burner  and  pool  fires  are  rather  simpler,  and  will  not  be  described  here. 

3.2  The  Flame 

The  flame  is  modeled  by  a homogeneous,  isothermal  gray  cone  of  gas,  of 
semi-apex  angle,  ij*.  When  the  fire  is  starved  for  oxygen,  or  when  near  fuel 
exhaustion,  the  cone  becomes  squatter,  thereby  emitting  less  radiation.  In 

Recalculation  in  1985  showed  that  y is  closer  to  the  stoichiometric  value, 
9.85. 
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When  oxygen 


the  absence  of  other  limiting  factors  the  cone  angle  is  \p 


30°. 


o 


starvation  or  burnout  reduce  the  burning  rate,  then 


(8) 


The  flame  gas  temperature  is  taken  to  be  1260°K,  and  its  absorption 


[15]  and  Orloff  [16].  Radiation  per  unit  area  to  the  flame  base  is  then  given 
by 


where  V = cone  volume,  A = bounding  area,  and  £ = shape  factor  (or  order 
0.9).  Radiation  from  the  flame  to  other  targets  is  complicated  by  the  fact 
that  the  slab  on  which  it  "sits"  may  shade  targets  lying  below  it. 
Expressions  for  the  flux  to  ceiling,  targets  (horizontal  or  vertical)  and  that 
absorbed  by  the  hot  layer,  are  given  in  references  [3]  and  [13] . 


The  plume  model  used  is  that  due  to  Morton,  Taylor  and  Turner  [17]  with 
the  heat  source  a virtual  point  below  the  surface  of  the  fuel  at  a distance 
such  that  the  plume  radius  at  the  surface  matches  the  fire  radius: 


coefficient  k = 1.55  m in  agreement  with  measurements  made  by  Markstein 


(7) 


3.3  The  Fire  Plume 


x = R/l . 2a 


(9) 
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Here,  a “ the  entrainment  coefficient,  taken  to  be  0.1 
the  plume  into  the  hot  layer  is 


The  mass  flow  from 


m 

p 


where 


p = ambient  air  density  (1.177  kg/m^  at  300°K) 

3 

b = 1.2  ahp  is  the  plume  radius  at  the  hot-cold  layer  interface 

hp  = hr  - hf  - hL  + x is  the  height  of  the  plume  above  the 

virtual  source.  Note  that  the  plume  height  shown  in  figure  1 

is  h = h'  - x. 

P p 

hr  = room  height 

h^  = height  of  the  fire  base 

hL  = depth  of  the  hot  gas  layer 
1/3 

u = (CQ/hp)  is  the  axial  flow  velocity  at  the  hot-cold 

layer  interface 
1/3 

U£  = (Cq/x)  is  the  axial  flow  velocity  at  the  fire  base 
C0  = 25  gQ/48  *o2CpTaPa 
Ta  = ambient  temperature 
Cp  = ambient  air  specific  heat 
g = acceleration  of  gravity 
Q = rate  of  heat  release  in  the  plume 
a = the  entrainment  coefficient  =0.1 


Enthalpy  is  carried  into  the  hot  layer  by  the  plume  at  the  rate 


E = m C T + Q-  E (16) 

P P P a pr 
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where  E is  t 
pr 

percent  of  Q). 


is  the  energy  lost  by  the  plume  via  radiation  (of  the  order  of  35 


3.4  The  Hot  Layer 


In  common  with  most  zone  models,  two  gas  layers  are  assumed.  The  upper 
gas  layer  is  homogeneous  and  isothermal,  and  the  lower  layer  consists  of 
ambient  air.  The  interface  between  the  two  layers  is  assumed  to  be  a hori- 
zontal plane.  In  reality,  of  course,  there  is  seldom  a distinct  homogeneous 
"upper  layer" — the  temperature  and  species  concentrations  vary  from  point  to 
point;  a first  simplification  of  reality  would  be  to  assume  a single  smoothly 
varying  vertical  profile  for  all  locations  within  the  room  (except  within  the 
fire  plumes).  The  step-function  profile  assumed  for  current  two-layer  models 
is  a further  simplification,  but  it  seems  adequate,  so  far. 

An  energy  balance  for  the  hot  layer  yields 


(17) 


where  E = the  energy  added  from  the  fire  plumes 


EE  = the  energy  carried  out  through  the  vents  by  mass  flow 
v ^ 

E^  = the  energy  convected  into  the  walls 


ET  „ = the  net  loss  of  energy  by  radiation 
L K. 


Similarly,  the  layer  mass  changes  according  to 


• _ • • 

til  = £ m — m 


(18) 


L 

P 
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where 


nip  = the  mass  added  from  the  fire  plumes 

= the  rate  at  which  mass  leaves  the  room  via  the  vent(s) 

The  layer  mass  and  energy  are  then  found  by  integrating  equations  [17]  and 
[18]. 


Since  we  can  also  write 


EL  = Wl  = KV  Cl\  ■ »LAhLCLTL  (19) 

(where  A is  the  cross-sectional  area  of  the  room,  is  the  specific  heat  of 

the  hot  layer  gases,  h^  is  the  layer  depth,  and  its  (mass)  density),  we 

have 


hL  = 


el/aVltl  = 


et  /act  p t 

L w a a 


(20) 


where  the  second  equality  follows  from  the  ideal  gas  law.  If  we  ignore  the 
(small)  difference  between  and  Cp  (the  specific  heat  of  ambient  air,  at 
constant  pressure),  it  is  clear  that  the  layer  depth  is  a function  only  of  the 
layer  energy.  It  also  becomes  clear  that  the  total  gas  energy  in  the  room  is 
(essentially)  constant:  as  energy  is  added,  enough  gas  is  pushed  out  by  gas 
expansion,  to  leave  the  total  energy  unchanged.  Finally,  the  layer 
temperature  is  found  from 


tl  ■ 


el/clml 


(21) 
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The  oxygen  concentration  in  the  layer  is  readily  found  as  the  solution  of  the 
equation 


m = 0.2318  T dl  - y 1 - Y m (22) 

ox  L L s dj  o u 

where  Y = m /nr  (23) 

o ox  L 

is  the  mass  fraction  of  oxygen  in  the  hot  layer  and  y is  the  stoichiometric 

s 

air/fuel  mass  ratio  (9.85  for  polyurethane);  0.2318  is,  of  course,  the  mass 
fraction  of  oxygen  in  air,  and  m^ , m . and  m^  have  all  been  defined  already. 
As  for  the  other  chemical  species,  each  burning  material  is  assumed  to 
generate  a constant  fraction  of  CO,  CO2,  ^0,  soot,  and  hydrocarbons  per  unit 
mass  of  fuel  pyrolyzed.  Values  for  these  fractions  are  taken  from  measure- 
ments [18]  and  the  species  concentrations  in  the  hot  layer  thus  become 
calculable.  These  concentrations,  in  turn,  permit  evaluation  of  the  layer's 
absorptivity,  k,  which  is  then  used  in  various  radiation  calculations.  Two 
alternative  routines  are  used  to  evaluate  k.  The  simpler  one  assumes  that  all 
the  absorption  is  due  to  (grey)  soot  absorption: 


k = 265  Y (24) 

s 

where  Y = m /mT  (25) 

s s L 

is  the  soot  mass  fraction.  The  more  complicated  one  also  takes  into  account 
the  band  absorption  of  CO2  and  l^O. 
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3.5  The  Vents 


The  mass  flows  through  each  vent  are  calculated  according  to  a hydraulic 
model  [19]  , with  an  efflux  coefficient  of  0.68.  The  pressure  at  the  floor  in 
the  center  of  the  room  is  used  as  an  auxiliary  variable.  As  a result  of 
assuming  the  gas  layers  to  be,  separately,  homogeneous,  the  pressure  differ- 
ence across  a vent  is  a piecewise  linear  function  of  height.  The  vent  is 
divided  into  horizontal  strips  corresponding  to  each  linear  variation,  and 
include  any  place  where  Ap  drops  to  zero  (a  neutral  axis)  as  a strip  boundary. 

A.  U 

Then  for  the  icn  strip,  the  pressure  drop  at  the  top  of  the  strip: 


(26) 


where  Ap.  is  the  pressure  drop  at  the  bottom  of  the  strip,  p is  the  ambient 

1 3 

air  density,  p that  of  the  hot  layer,  and  Ah^  the  height  of  the  strip.  The 
mass  flow  through  that  strip  then  is 


where 


m. 


„ b + /ab  + a 

“ i = - 

/ b + /a 

5 | AP± | , b = |Ap1+1 


(27) 

(28) 


and  Gi  = (sign  Ap^  (2/3)  CdB  At^  /2  gp^p 


(29) 


Here  C . = 0.68  and  B is  the  width  of  the  vent.  When  Apc1  is  correctly 
a r loor 

chosen,  all  the  flows  are  such  that  mass  balance  for  the  room  is  satisfied: 


m 

out 


m 


in 


h 

C T 
P a 


- E 
f 


ra„ 


(30) 
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3.6  Radiation 


There  are  eight  radiation  subroutines  which  calculate  the  most  important 
radiation  exchanges  among  the  flame(s),  objects,  hot  layer,  and  hot 
wall/ceiling.  Some  are  quite  simple,  like  equation  (7).  Most  others  are 

complex,  and  will  not  be  listed  here;  the  interested  reader  should  see 
reference  [3]. 


The  convective  heat  exchanges  between  the  hot  layer  and  walls,  and  that 
between  layer(s)  and  objects,  are  found  assuming  a simple  (but  temperature- 
dependent)  film  coefficient.  For  heat  transfer  to  ambient  air,  Mark  5 takes 


and  Tg  is  the  temperature  of  the  surface  being  cooled  (outside  of  walls,  or 
objects  in  the  room).  For  heat  transfer  from  the  hot  layer, 


3.7  Convective  Heat  Transfer 


(31) 


where 


(32) 


where 


h . = min 

l 


[50,5  + 0.45  (Tl  - Ta)] 


(33) 


In  the  steady  state,  the  heat  loss  rate  to  the  walls  and  ceiling  therefore  is 


(34) 


®LW  ~ \ \ 


where  = 2 (L  + W)  hL  + LW  - Av 


(35) 


where  L and  W are  the  length  and  width  of  the  room,  and  Ay  is  that  part  of  the 
vent  area  intersected  by  the  hot  layer.  Here  T0  of  equation  (32)  is  taken  to 
be  T^,  the  inner  wall/ceiling  temperature,  assumed  to  be  uniform.  As  the  fire 
grows,  however,  the  layer  thickens  and,  in  order  to  (artificially)  maintain  a 
uniform  temperature,  extra  heat  has  to  be  supplied  to  the  newly  exposed  wall 
areas.  This  must  be  done  such  that  the  newly  heated  wall  immediately  acquires 
the  same  temperature  profile  throughout  its  thickness  as  the  rest  of  the 
heated  wall.  Hence  for  the  growing  layer,  equation  (34)  is  replaced  by 

- Vw  + 2 <L  ♦ w>  t"  \ (\  > 0)  (30) 

where  q"  = / (q”  + q^w)  dt  (37) 

o 

is  the  energy  per  unit  area  stored  in  the  wall,  and  q"  is  given  by  equation 

AW 

(31),  with  T = T . , £ .It  follows  from  equation  (20)  that 

s outside  wall  surface 

hL  = hLEL/EL  (38) 

When  the  layer  decreases  in  depth,  we  evidently  must  use  equation  (34)  r.ither 
than  equation  (36). 
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3.8  Conductive  Heat  Transfer 


Walls  and  objects  are  assumed  to  be  finite-thickness  slabs  insulated  at 
their  edges.  Thus  the  one-dimensional  heat-conduction  equation  can  be  used. 


3T 

3T 


(39) 


is  solved  numerically  with  a simple  algorithm.  Although  this  finds  the 
temperature  profile  through  the  object  as  a function  of  time,  most  of  that 
information  is  neglected,  as  we  are  (currently)  only  interested  in  the  surface 
temperatures  and  the  integral  of  the  profile. 


This  completes  the  brief  outline  of  the  physics  included  in  the  model. 


The  numerical  methods  by  which  the  resulting  equations  are  solved  will 
not  be  discussed  since  this  is  done  in  detail  in  CFC  III  [9]  and  CFC  V [13]. 
For  those  readers  who  are  not  familiar  with  these  papers,  we  might  mention 
that  two  alternative  algorithms  are  used:  a Gauss-Seidel  successive  substi- 

tution method  is  normally  used;  when  convergence  is  not  achieved  with  that 
method,  a multivariate  Newton  technique  is  used.  This  combination  works  very 
well,  95  percent  of  the  time. 


4 . THE  TESTS 


In  1977,  Factory  Mutual  Research  Corporation  carried  out  [20]  a series  of 
eight  well-instrumented  full-scale  room  fire  tests,  against  which  the  Mark  5 
program  can  be  checked.  The  eight  tests  were  systematic  variations  of  a 
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a standard  room.  The  standard  configuration  is  shown  in  figure  2;  it  is  a 
2.44  by  3.66  by  2.44  m room  (8*  x 12'  x 8')  with  an  open  doorway  0.76  x 2.03  m 
(30"  wide  by  80"  high);  a slab  of  polyurethane  (P.U.)  foam,  measuring  1.52  x 
1.52  x 0.10  m (5'  x 5'  x 4")  is  in  one  corner;  its  top  surface  is  61  cm  from 
the  floor  (this  models  a foam  mattress).  A P.U.  foam  target  1.22  x 0.30  x 
0.10  m (4'  x 1*  x 4")  is  in  the  near  adjacent  corner;  its  top  surface  is  25  cm 
(10")  above  that  of  the  primary  fuel  slab.  The  primary  slab  is  ignited  at  its 
center,  and  the  fire  followed. 

Tests  #0,  1,  and  2 were  all  meant  to  be  identical — the  standard 
configuration.  Test  3 had  the  doorway  narrowed  to  19  cm;  test  4 added  a 
window;  test  5 had  a pool  of  PMMA  beads  replacing  the  foam;  test  6 had  a 
window  instead  of  a door,  and  test  7 had  the  door  reaching  up  to  the  ceiling. 
Comparing  the  results  of  calculations  with  the  results  of  these  experiments 
shows  that  the  "predictions"  vary  from  excellent  to  fair  for  most  of  the 
variables.  Although  comparisons  have  been  made  for  all  the  tests,  for  brevity 
only  two  of  these  will  be  discussed:  the  standard  configuration  and  test 

#6.  These  show  the  best  and  poorest  agreements  between  the  predictions  of  the 
model  and  the  actual  fires. 


5.  RESULTS 

The  rate  of  mass  loss  of  the  fuel  (while  it  is  burning)  is  shown  for  runs 
0,  1,  and  2,  in  figure  3.  These  are  curves  A,  B,  and  C,  respectively.  We  see 
that  A and  B agree  very  well,  although  there  is  a good  deal  of  statistical 
fluctuation  for  t < 150  sec.  For  test  #2  (on  July  8)  fire  development  was 
more  rapid,  as  can  be  seen  from  curve  C — it  looks  as  if  ignition  occurred 
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about  40  seconds  early.  Note  that  these  are  semilog  plots,  as  is  also  true 
for  figs.  6,  11,  12,  and  13. 

The  thickness  (or  depth)  of  the  hot  gas  layer  for  each  of  these  first 
three  tests  is  shown  in  figure  4.  Up  to  t = 200  sec,  curves  A and  C (tests  0 

and  2)  track  together  very  well,  whereas  curve  B (the  dotted  curve, 

corresponding  to  test  1)  begins  its  descent  some  50  sec  earlier;  this 
contrasts  with  the  behavior  in  fig.  3,  where  ignition  appeared  to  have 
occurred  earlier  for  test  2,  rather  than  test  1.  Beyond  200  sec,  B and  C 
switch  - that  is,  curve  C leads  the  other  two  by  about  40  sec  (consistent  with 
the  mass-loss  rate  shown  in  fig.  3),  while  A and  B track  together  quite  well. 
This  holds  up  to  t = 310  sec;  after  that,  there  are  substantial  deviations, 
curve  B displaying  marked  oscillations  whose  nature  is  not  clear.  These  two 
figures  serve  to  show  how  much  variation  there  can  be  in  reasonably  well- 
controlled  experiments — so  that  we  should  not  expect  calculated  results  to 
agree  any  better  than  this,  with  any  one  experimental  fire. 

Test  #1  (July  6)  was  used  as  the  standard  for  comparison,  partly  because 

a roof  leak  appeared  during  test  #0,  and  it  is,  therefore,  suspect.  The 

target  ignited  310  seconds  after  ignition  of  the  primary  slab,  and  the 
sprinkler  was  turned  on  at  t = 326  sec. 

Figure  5 shows  how  differences  in  the  handling  of  the  data  can  also  lead 
to  rather  different  results.  The  solid  curve  shows  "the  temperature"  T^(t)  of 
the  hot  layer  as  a function  of  time,  as  obtained  by  taking  an  average  of  the 
readings  of  8 thermocouples  located  at  different  places  in  the  room,  all  at  a 
height  1.71  m from  the  floor — (70  percent  of  the  distance  to  the  ceiling). 
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This  was  the  value  of  T^Ct)  used  by  Factory  Mutual  in  their  data  analysis.  We 
made  a different  evaluation  of  T^Ct),  using  the  temperature  profiles  as  a 
function  of  height  at  three  places  in  the  room;  this  gave  the  dashed  curve. 
As  can  be  seen,  the  latter  method  gives  a peak  temperature  150°C  lower  than 
the  FMRC  method  does. 

6.  PREDICTION  VS  EXPERIMENT 


Figure  6 shows  the  pyrolysis  rate  for  test  1H  (July  6).  The  calculation 
generally  agrees  with  experiment,  though  the  maximum  is  some  50  percent  high. 
Figure  7 shows  the  hot  layer  temperature  for  the  same  test.  As  we  see,  the 
calculation  agrees  quite  well  with  experiment,  out  to  310  sec.,  but  then  falls 
below  the  experimental  value.  Incidentally,  no  meaningful  comparisons  can  be 
made  after  t = 326,  when  the  sprinkler  was  turned  on — the  program  does  not  yet 
simulate  extinguishment.  Note  that  the  calculated  peak  temperature  is  low 
even  though  the  burning  rate  is  over- predicted  at  this  point. 

In  figure  8 the  layer  height  calculation  is  shown  to  agree  very  well  with 
experiment  out  to  t = 290;  the  subsequent  disagreement  occurs  in  part  because 
in  Mark  5,  the  layer  cannot  descend  below  the  level  of  the  base  of  the  fire. 
It  is  possible,  moreover,  that  late  in  the  fire  the  interface  between  the  hot 
and  cold  layers  is  not,  in  fact,  a plane:  it  may  be  that  the  surface  rises 

locally,  near  the  fire,  and  is  lower  elsewhere;  because  of  the  way  in  which  it 
was  obtained  the  data  would  only  reflect  the  lower  part  of  the  surface. 

The  outflow  rate  of  hot  gas  is  shown  in  figure  9 — again,  the  calculation 
is  in  excellent  agreement  with  experiment  out  to  t = 300  sec.;  beyond  that, 
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the  hotter  layer  produces  greater  buoyancy  (recall  figure  7),  and  therefore 
greater  flow.  The  target  temperature  (not  shown)  also  agrees  very  well  with 
calculation. 

Figure  10  shows  that  the  oxygen  concentration  in  the  layer  is  reasonably 
well  calculated.  The  fact  that  the  real  minimum  is  lower  than  the  calculated 
minimum  is  of  course  consistent  with  the  temperature  and  mass  flow  curves, 
though  a possible  alternative  explanation  will  be  given  when  we  discuss  test 

#6. 


Other  species  concentrations  are  shown  in  figure  11;  the  calculated  CO2 
concentration  agrees  very  well  with  experiment,  but  the  CO  does  not. 

The  radius  of  the  fire  as  a function  of  time  is  well  predicted  by  the 
program,  as  is  the  radiative  flux  at  the  center  of  the  floor  (most  of  the 
radiation  comes  from  the  hot  layer);  see  figure  12.  (The  calculated  flux 
rises  too  slowly  beyond  t = 305  because  the  calculated  layer  temperature  rises 
too  slowly  there.) 

In  sum,  the  standard  run  reproduces  the  results  of  tests  0 and  1 quite 
well,  the  principal  disagreement  being  the  failure  to  reproduce  the 
temperature  peak,  and  the  consequences  which  flow  therefrom. 

Next,  consider  test  #6,  where  the  doorway  was  replaced  by  a window  about 
1 m high,  m^,  T^,  h^,  m^,  Tw  and  y(02)  will  be  used  for  comparison. 
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The  pyrolyzation  rate,  mf,  (see  fig.  13)  is  in  excellent  agreement  up  to 
t = 300,  but  the  calculated  peak  again  continues  to  rise  beyond  the  experi- 
mental peak.  In  spite  of  that,  the  calculated  temperature,  T^,  rises  more 
8 lowly  than  the  experimental  value,  as  seen  in  figure  14,  and  dies  out  much 
faster,  beyond  t = 360. 

The  layer  depth  in  the  room,  hL,  is  shown  in  figure  15,  and  is  seen  to 
agree  reasonably  well  with  experiment  between  t = 90  and  t = 230.  Beyond  that 
point,  the  layer  is  again  quite  a bit  deeper  than  calculated.  The  disagree- 
ment is  bad  after  t = 360,  just  as  it  is  for  the  temperature;  however,  the 
experimental  values  are  themselves  suspect  out  here. 

The  resulting  gas  outflow  rate,  m^,  is  shown  in  figure  16,  and  is  in 
surprisingly  good  agreement  with  experiment  out  to  t = 210.  The  first  calcu- 
lated maximum,  at  280  seconds,  is  due  to  oxygen  starvation:  the  layer  is  so 
low  that  not  enough  fresh  air  is  drawn  in  by  the  plume  to  permit  burning  as  in 
the  open.  The  second  maximum  (at  t = 350)  is  due  to  the  beginning  of  burnout 
of  the  primary  fuel.  The  third  peak,  at  t = 380,  is  not  understood. 

In  figure  17,  we  show  the  average  temperature  of  the  ceiling  and  upper 
part  of  the  walls  (that  part  in  contact  with  the  hot  layer)  Ty.  As  could  be 
expected,  the  curves  are  similar  to  the  layer  temperature  curves  (see  fig. 
14) — that  is,  the  calculated  wall  temperature  lags  behind  the  measured  temper- 
ature by  about  15  sec,  out  to  t = 290,  then  falls  behind  more  dramatically. 

The  oxygen  concentration  at  the  window,  y(02),  is  shown  in  figure  18  as  a 
function  of  time  for  probes  at  different  heights.  Note  that  in  this  figure 
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♦-he  heights  refer  to  distance  down  from  the  ceiling.  Thus  the  curve  marked 
0.2H  gives  the  oxygen  concentration  20  cm  below  the  ceiling.  The  stratifica- 
tion of  the  layer,  which  was  not  evident  in  test  #1,  is  brought  out  most 
strikingly,  here.  It  is  clear  from  examination  of  the  test  data  that  there 
may  well  have  been  variations  along  the  length  of  the  room.  (The 
stratification,  incidentally,  has  important  implications  for  burning  in  the 
vitiated  layer.)  Were  the  layer  to  be  well  stirred,  the  mean  oxygen 
concentration  would  evidently  be  a weighted  average  of  the  experimental 
values.  This  average  was  not  taken,  but  it  is  probable  that  the  calculated  O2 
concentration  would  agree  reasonably  well  with  such  a mean. 

Finally,  the  fact  that  the  oxygen  concentration  at  the  fire  base  is 
ambient  throughout  most  of  the  burn  in  spite  of  the  fact  that  the  layer  has 
fallen  below  the  base  of  the  fire,  is  not  consistent  with  the  assumption  that 
the  layer  is  flat  everywhere.  It  is  much  more  likely  that  it  rises  locally, 
near  the  fire,  as  mentioned  earlier. 

7.  DISCUSSION 


The  Mark  5 growing  fire  algorithm  predicts  the  radius  of  the  fires  quite 
well.  The  same  is  true  of  the  pyrolysis  rates  early  in  the  fire;  however,  the 
predicted  peaks  tend  to  be  too  high,  and  the  resulting  burnout  too  swift. 

The  vent  routine  is  quite  good,  and  gives  mass  flow  rates  in  excellent 
agreement  with  measurement  provided  the  layer  temperature  is  well  predicted. 
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The  plume  model  we  use,  though  crude,  seems  satisfactory,  judging  only  by 
the  amount  of  air  entrained  and  the  resulting  oxygen  and  CO2  concentrations  in 
the  layer.  The  flux  from  the  hot  layer  to  the  floor  agrees  well  with  experi- 
ment, suggesting  that  the  layer  emissivity  is  reasonably  well  estimated,  and 
therefore  that  the  smoke  concentration  is  reasonably  well  predicted. 

On  the  other  hand,  the  wall  temperatures  are  a bit  low,  suggesting  that 
radiation  from  the  flame  plus  convection  from  the  layer  are  under-estimated. 
This  suggests  that  the  flame  model  is  not  quite  adequate.  This  impression  is 
reinforced  when  we  note  that  the  flux  to  the  targets  is  low  in  some  cases  (but 
not  in  others). 

The  burnout  of  the  fuel  was  made  _ad  hoc , and  turns  out  to  be  poor.  This 
is  not  of  central  significance  at  the  moment,  but  clearly  has  to  be  improved. 

Because  CO  is  a major  hazard  of  fire  the  failure  of  the  model  to  predict 
its  concentration  well  is  a major  weakness.  CO  production  in  diffusion  flames 
is  not  well  understood  and  the  model  reflects  this.  Fortunately,  research 
currently  in  progress  shows  promise  that  a better  understanding  of  the  physics 
and  chemistry  of  CO  production  will  soon  be  available. 

Perhaps  the  most  puzzling  question  is  why  the  layer  temperature  is  well 
predicted  for  most  of  the  burn,  but  underpredicted  late  in  the  fire.  There 
are  at  least  five  possible  reasons  for  this: 
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1.  Mark  5 may  overestimate  the  rate  of  heat  loss  by  the  layer. 

2.  There  may  be  burning  in  the  layer  itself. 

3.  The  plume  model  may  be  inadequate,  calculating  more  ambient  air 
entrainment  than  is  found  in  practice. 

4.  Mark  5 neglects  heating  of  the  lower  layer;  this  may  be  important. 

5.  Air  may  be  drawn  up  around  the  base  of  the  flames.  Mark  5 may  limit 
burning  due  to  oxygen  starvation  prematurely. 

Let  us  consider  these  in  order:  //I  is  very  unlikely  since  it  would  lead  to 

wall  temperatures  higher  than  observed,  rather  than  lower.  Also,  energy  loss 
through  the  vent  is  already  underestimated.  #2  is  unlikely,  as  no  thermo- 
couple gave  a reading  anywhere  near  flame  temperature.  (However,  the  thermo- 
couple response  times  may  be  too  long.)  #3  is  possible,  though  not  very 
likely.  Heating  of  the  lower  layer,  however,  ( #4 ) certainly  occurs — both 
because  the  floor  is  heated  by  radiation,  thereby  heating  the  air,  and  by 
entrainment  of  the  hot  gases  by  the  ambient  air  coming  in  through  the  vents. 
This  will  make  the  upper  layer  hotter,  too. 

The  other  effect  which  has  been  neglected,  is  the  raising  of  the  hot/cold 
interface  at  the  fire.  This  would  delay  oxygen  starvation,  and  allow  more 
vigorous  burning  to  take  place.  These  are  matters  to  be  investigated,  both 
theoretically  and  experimentally. 
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8.  CONCLUSIONS 


The  overall  performance  of  the  Mark  5 simulations  is  satisfactory,  both 
in  its  numerical  stability  and  convergence,  and  in  the  predictions  it  makes. 
There  are  a small  number  of  inadequacies,  which  are  being  investigated;  there 
is  every  reason  to  believe  these  difficulties  will  be  overcome. 


In  spite  of  all  their  complexities,  therefore,  it  appears  that  fires  can 
be  successfully  simulated  by  mathematical  models,  and  it  is  clear  that 
continued  efforts  will  lead  to  continued  progress. 
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